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I. INTRODUCTION 

One of the principal deficiencies in our understanding 
of liquid mixtures is the nature of the link between the 
microscopic description of the system and its macroscopic 
phase behaviour. For simple single component fluids, the 
phase diagram topology is relatively insensitive to the mi- 
croscopic properties of the molecules and exhibits (even 
in systems with strong anisotropics or long ranged inter- 
actions) a liquid-vapour first order line terminating at a 
critical point. By contrast, in binary mixtures, the in- 
terplay between the constituent components leads to a 
wealth of intriguing phase behaviour depending on 
the relative sizes of the molecules and the strengths of 
their interactions. 

Although the various possible phase diagram topolo- 
gies have been placed into a number of categories [Q, 
it is not well understood (even at the mean field level) 
precisely which microscopic features are responsible for 
yielding a given topology. Also unclear is the extent 
to which critical fluctuations affect the structure of the 
phase diagram i.e. whether the neglect of correlations 
in many analytical theories yields qualitatively (as well 
as quantitatively) incorrect phase diagrams. The task of 
accurately and reliably deriving the full phase behaviour 
of a fluid mixture from knowledge of its microscopic in- 
teractions therefore remains a great challenge. 

Evidently, an accurate description of the phase be- 
haviour of a simple binary fluid model, would provide 
a useful benchmark against which current and future 
liquid-state theories could be tested. In the present work, 
we furnish such a description by means of Monte Carlo 



simulations of a simple continuum model, the results of 
which we compare with mean field calculations. For rea- 
sons of computational tractability, we consider a sym- 
metrical binary fluid model, i.e. one in which the two 
pure components A and B are identical and only the in- 
teractions between particles of dissimilar species differ. 
Notwithstanding its simplicity, however, the model turns 
out to reveal a rich variety of interesting phase behaviour. 



II. BACKGROUND 

The phase diagram of a symmetrical binary fluid mix- 
ture is spanned by three thermodynamic flelds (T, /x, h) , 
where T is the temperature, fi is the overall chemical po- 
tential coupling to the total density, and h is an ordering 
field coupling to the relative concentrations of the two 
fiuid components which we assume are allowed to fluctu- 
ate. In this work we shall restrict our attention to the 
phase behaviour in the symmetry plane h — 0, i.e. we 
stipulate that on average the numbers of A and B par- 
ticles are equal. Additionally we shall assume that simi- 
lar species interactions are energetically more favourable 
than dissimilar species interactions. This latter condition 
provides for a consolute point (critical demixing transi- 
tion) at some finite temperature Tc. For temperatures 
T < Tc, there is coexistence between an A- rich liquid 
and a i?-rich liquid, while for T > the system com- 
prises a homogeneous mix of A and B particles. Precisely 
at Tc, the system will be characterised by strong critical 
concentration fluctuations 



between the ^-rich and B-rich phases. Such a demixing 
transition is analogous to that occurring at the critical 
point of a simple spin-^ Ising model. The difference for 
an off-lattice fluid, however, is that the demixing temper- 
ature depends on the density. Consequently one obtains 
a critical line of consolute points Tc{p) [or Tc(/i)], which 
is commonly referred to as the 'A line'. 

In addition to exhibiting consolute critical behaviour, 
binary fluids can also exhibit liquid-vapour (LV) coexis- 
tence, in much the same way as does a single component 
fluid. It transpires, however, that the LV phase behaviour 
of binary mixtures is considerable richer than that of sim- 
ple fluids. This difference is traceable to the additional 
ingredient of concentration fluctuations, which couple to 
the density fluctuations and which can radically alter the 
LV phase behaviour. Since concentration fluctuations are 
strongest on the A line, one expects that alterations to 
the LV coexistence behaviour will be greatest where this 
line approaches the LV coexistence curve. 

Perhaps not surprisingly, binary fluids mixtures are not 
the only fluid systems in which first order phase coexis- 
tence behaviour is influenced by the proximity of a criti- 
cal line. The earliest sightings of such effects appears to 
have been in analytical studies of various lattice-based 
fluid models |^-|^. Some time later, a detailed Landau 
theory study of a model for sponge phases in surfactant 
solution 1^, revealed a rich variety of first order phase 
behaviour as the path of the A line was varied. More 
recently, similar behaviour was uncovered in extensive 
mean field and density functional theory investigations of 
a number of symmetrical continuum fluid models, namely 
the classical Heisenberg spin fluid [f7|-|lO|, a dipolar fluid 
model and the van der Waals-Potts fluid 

Despite dealing with ostensibly quite distinct mod- 
els, the gross features of the mean field phase be- 
haviour emerging from these studies appears to be essen- 
tially model- independent. This behaviour is illustrated 
schematically in fig. |l| and involves three possible LV 
phase diagram topologies, depending on the path of the 
critical line relative to the LV line. To describe this be- 
haviour we shall employ the language of the symmetrical 
binaryfluid. In so doing, we anticipate the result of sec- 
tions [II and namely that the same scenario is played 
out in this case too. Of course, to obtain the correspond- 
ing behaviour for other systems, eg. the magnetic or 
dipolar fluids, one need only substitute the appropriate 
nomenclature eg. 'mixed fluid' — > 'paramagnetic fluid'. 
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FIG. 1. Schematic representation of the three types of 
phase diagram for a symmetrical binary fluid mixture in the 
density-temperature plane, as described in the text. The full 
curve is the first order liquid- vapour coexistence envelope, 
while the dashed curve is the A line of critical demixing tran- 
sitions. 



Fig. |l|a depicts schematically the mean field phase dia- 
gram obtained when the model parameters are chosen 
such that the A line approaches the first order phase 
boundary well below the liquid-vapour critical point. In 
such a situation, the A line intersects the LV line at a 
critical end point (CEP). At the CEP, a critical liquid 
coexists with a non-critical vapour. Below the CEP tem- 
perature one finds a triple line in which a vapour coexists 
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with an A-rich liquid and a _B-rich liquid. Owing to the 
symmetry, these two liquids have the same density. 

Alternatively, for different model parameters, the A 
line may intersect the LV line at the liquid- vapour critical 
point [fig|^c]. Under such conditions, phase coexistence 
between the vapour and the mixed fluid is preempted by 
the demixed fluid phase. One then obtains a tricritical 
point |]l^ , |l5| in which three phases (a vapour, an A-rich 
liquid, a i?-rich liquid) simultaneously become critical. 

The intermediate situation is shown in fig. P?, and oc- 
curs when the A line approaches the LV line at a tem- 
perature somewhat (but not greatly) below the liquid- 
vapour critical temperature. In this case, the phase di- 
agram combines the features of the previous two cases. 
One finds a triple point in which a vapour coexists with 
a mixed liquid at intermediate density and a demixed 
liquid of higher density [p^ . Above the triple point tem- 
perature, a demixed vapour and a demixed liquid coexist 
at low and moderate densities, becoming identical above 
the liquid- vapour critical point. At higher densities, a 
mixed liquid and the demixed liquid coexist, becoming 
identical at a tricritical point. 

That the scenario described above is generic to a range 
of apparently distinct fluid models (eg. dipolar, magnetic 
and binary fluids), is perhaps slightly surprising at first 
sight. On closer examination, however, it becomes clear 
that the model differences are only skin-deep. All the 
systems in which this behaviour has yet been identified, 
can, in essence, be regarded as fluids in which each par- 
ticle carries an internal degree of freedom, eg. a spin 
or dipolar moment. The symmetrical binary fluid model 
shares this behaviour because the particle species label 
is analogous to a two-state 'spin' variable. 

Notwithstanding the substantial body of analytical ev- 
idence supporting the scenario of phase behaviour shown 
in flg. |l], it must necessarily be regarded as somewhat 
tentative given the notorious inability of mean fleld the- 
ories to account accurately for critical behaviour below 
the upper critical dimension. In view of this, independent 
corroboration by computer simulation is clearly desirable 
and necessary. In recent times, such studies have indeed 
started to appear. 

The first simulations to study the confluence of a crit- 
ical demixing line and a first order LV line sought to 
elucidate the behaviour when the intersection occurs at 
the LV critical point, i.e. at a tricritical point JcL fig. |l|c]. 
Investigations of a two-dimensional spin fluid |17[ demon- 
strated that the tricritical point properties are identical 
to those of the 2D Blumc-Capel model, as one might ex- 
pect on universality grounds. Other studies, this time 
focusing on the 3D classical Heisenberg fluid (l^Jl^ and 
spin-i quantum fluids ||l9|-|2l[|, mapped the LV phase en- 
velope and the A line around the tricritical point and 
compared the results with mean fleld calculations. Mod- 
est agreement was found. 

As regards the situation shown in flg. ^ and there 
is a still a paucity of simulation data. This is presum- 
able traceable to the practical difficulties associated with 



studying first order phase coexistence deep within the 
two phase region. The basic problem is the ergodic (free 
energy) barrier to sampling both coexisting phases in a 
single simulation. This barrier arises because the phase 
space path leading from one pure phase to another nec- 
essarily passes through interfacial configurations of large 
free energy, having a concomitantly low statistical weight. 
Such configurations occur only very rarely in a standard 
Monte Carlo (MC) simulation. Fortunately, however, a 
recently introduced biased-sampling technique known as 
multicanonical preweighting |22| allows one to negotiate 
this barrier, and thus obtain accurate estimates of co- 
existence properties The efficacy of the method 
for investigating fluid phase coexistence was flrst demon- 
strated in [|4j. Very recently, it has also been employed 
to study critical end point behaviour in a Lennard- Jones 
binary fluid model [cf. flg. |^a]. In this study, the 
intersection with the A line was shown to engender a sin- 
gularity in the first order phase boundary-in accord with 
earlier theoretical predictions j2^ . For the LV phase en- 
velope, this singularity is manifest as a bulge in the liquid 
branch density (as indicated schematically in fig. |^a). 

Thus while fragments of the picture of LV behaviour in 
symmetrical fiuids have been set in place by simulation, 
much clearly remains to be done before a reliable and 
comprehensive overview emerges. In particular, no evi- 
dence has yet been reported for the existence of the triple 
point behaviour shown schematically in fig. ^3. There has 
also been no systematic simulation study of the full range 
of LV phase behaviour for a single model. Consequently 
little reliable information exists concerning the manner 
in which one type of phase diagram evolves into another. 

In the present work we have attempted to remedy 
this situation by performing a systematic MC simulation 
study of the LV phase behaviour of a binary fluid mix- 
ture. The model studied has an interparticle potential of 
the square-well form: 
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U{r) 
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r < a 

a < r < 1.5a 
r > 1.5(7 



(2.1) 



Here r is the particle separation, J is the well depth 
or interaction strength, and a is the hard-core radius. 
In general there will be a number of different interac- 
tion strengths J^^,J^^,J^^ depending on the respec- 
tive species of the interacting particles. However, in the 
symmetrical case with which we shall be concerned, one 
has simply 



J^^{r) = J^^(r)- J(r) 
J^S(r) = SJir). 



(2.2) 



The parameter S = Jab/ J ^ 1 determines the degree 
to which interactions between dissimilar species are less 
favourable than similar species interactions. Since 5 is 
the only free model parameter, it controls the complete 
range of possible phase behaviour. 
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Using MC simulation we obtain the liquid-vapour 
phase behaviour of this system for a range of values of 
6. The results are compared with explicit mean field 
model calculations in order to assess the latter's ability 
to reproduce the actual phase behaviour. Landau theory 
calculations are also reported which furnish physical in- 
sight into the general mechanisms by which the coupling 
of density and concentration fluctuations engender the 
various types of observed phase behaviour. Additionally 
we calculate (both explicitely for our model and within 
Landau theory) the spinodal lines of the phase diagram. 
This reveals that for a certain range of S there exists a 
'hidden' or metastable binodal, the presence of which is 
expected to have implications for the system dynamics. 

The remainder of our paper is organised as follows. 
In section [11 we detail the explicit mean-field and gen- 
eral Landau theory calculations for the phase diagram 
as a fimction of 6. The Monte Carlo simulation results 
for the liquid-vapour phase behaviour are presented in 
section Finally, in section ^ we conclude by compar- 
ing the simulation and mean field results, and discussing 
prospects for future work. 



III. MEAN-FIELD CALCULATIONS 



this we introduce a two-state species variable Si taking 
the value Si = 1 (—1) when the ith particle is of type A 
[B). The total configurational energy can then be writ- 
ten 

(i<i) (i<f) 

+ (3-2) 

(i<j) 

To obtain the liquid-vapour phase envelope and the A 
line, we adopt an approach similar to that employed in 
references p^- |29|pT| , p^ , ^ . In the thermodynamic limit, 
the (Helmholtz) free energy density f{p, T) as a function 
of the number density p — N/V and temperature T, is 
given by: 

/(p,T) = hm ^ Intr [exp (-/3$^)] . (3.3) 

Now, within the mean field approximation, one assumes 
an interaction between an A-type particle and an effective 
field: 



m, 



(3.4) 



i>i 



In this section we present two complementary mean 
field studies of symmetrical binary mixtures. The first 
is an explicit investigation of the square-well binary mix- 
ture model, also studied by simulation in section |V|. The 
second is a general Landau theory treatment aimed at 
obtaining physical insight into the origin of the various 
phase diagram topologies. 



A. Explicit model calculations 

Let us consider a binary fluid in which the two particle 
species A and B interact via the symmetrical square- 
well potential of equations 2.1 and |2.2|. For analysis pur- 



poses it will prove useful to decompose this potential into 
a hard-sphere component plus a short ranged attractive 
part. To this end we rewrite equation 2.1 as 



U{r)=UHsir) + Jir), 



(3.1a) 



and an interaction between a B-type particle and an ef- 
fective field 



J+-{J-/N)J2s^^J+-J- 



(3.5) 



i>l 



Here, are effective potentials and m = {Na — 
Nb)/{Na + Nb)- Since the coordination number in the 
fluid is density dependent, we make the approximation: 



-P / a r g(r). 



(3.6) 



where the fluid correlation function g[r) is taken from 
the Percus-Yevick solution for hard spheres |^ . 
The mean field configurational energy is then 



1 ^ 

--J2[hA{l + s,) + hB{l-s,)] (3.7) 



where 



and 



Unsir) = oo 
UHsir) = 



J(r) = - J 
J(r) = 



r < a 
otherwise 



(7 < r < 1.5a 
otherwise. 



(3.1b) 



(3.1c) 



For a binary fluid, the interaction between two parti- 
cles depends on their respective species. To deal with 



from which the free energy at constant T follows as : 
fMFip) - lim \ntr [cxp {-P^^f)] 

V^OO jj V 



= fHs{p) - P — 



^™ ^ ■§ 1^ 2 cosh {pj^m) 
I p 



(3.8) 



Here fusip) is the free energy of a reference system com- 
prising a hard-sphere single component fluid. The third 
term of the right hand side of eqn. (3.8) is minimised for 
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m = tanh [/3J m] 



(3.9) 



As f(p) is not always a convex function of the density, we 
take the convex envelope in order to find the coexistence 
densities for the first order LV transition. 

We also wish to obtain the phase diagram in the fj,-T 
plane. To achieve this one needs to consider the grand 
potential f{p) + fip, with the chemical potential. Min- 
imising this yields the pressure: 

p{p, p) = nmv [-/(/>') + lip'] ■ (3.10) 

The coexistence chemical potential is then found by de- 
manding equality of both the chemical potential and the 
pressure in the coexisting phases. 

The resulting phase diagrams in the p-T plane are 
shown in fig. ^-|^, with the p-T phase diagrams shown 
as insets. For large values of 5 < 1 we find a LV coexis- 
tence region and a A line at high densities which intersects 
the LV line at a critical end point. The CEP induces an 
anomaly or kink in the liquid branch density which is 
clearly visible in fig. ||a. This anomaly is the mean field 
remnant of the specific heat-like singularity studied the- 
oretically and computationally in references |2^j2^. Be- 
cause the mean field specific heat exhibits a jump rather 
than a divergence at criticality, the anomaly takes the 
form of an abrupt change in the gradient of the liquid 
branch dpi/dT. If fiuctuations are taken into account, 
however, dpi/dT diverges at the CEP [^ . 

As 5 is reduced, the CEP anomaly grows until at 
around 5 = 0.7 a small peak emerges [cf. fig ^d]. The 
point at which this occurs con stitut es a tricritical end 



point as discussed in subsection |[II B| . On further reduc- 
tion of S, the peak develops until the situation shown 
in fig. ^ is attained. Here the first order coexistence 
envelope displays a triple point at which a vapour coex- 
ists with a mixed liquid at intermediate density and a 
demixed liquid of higher density. Above the triple point 
temperature, a vapour and a demixed liquid coexist at 
low and moderate densities, merging at the liquid- vapour 
critical point. At higher densities a mixed liquid and the 
two symmetrical demixed liquids coexist, becoming iden- 
tical at a tricritical point. 

If S is reduced further still, one reaches a point (for 
6 < 0.605), at which the triple point temperature equals 
the LV critical point temperature. Thereafter, the liquid- 
vapour critical point is lost and only a tricritical point 
remains, as shown for the case 6 = 0.57 in fig. |d. No 
further topological changes in the phase diagram are ob- 
served as S is reduced to zero. 




FIG. 2. Mean field phase diagrams in the p-T plane for 
various 5, as described in the text, (a) 5 = 0.72, (b) 5 — 0.70, 
(c) 5 = 0.65, (d) S = 0.57. The insets show the corresponding 
phase diagrams in the p-T plane. Full lines represent first 
order phase coexistence and dashed curves represent the A 
line. The demixed fiuid (DF), mixed fiuid (MF) and vapour 
(V) phases are also marked. 



Although the hquid-vapour critical point disappears 
from the equilibrium phase diagram for d < 0.605, it 
is interesting to note that it nevertheless remains in 
metastable form for some range of 6. This is illustrated 
in fig. ^ which shows the phase diagram for S = 0.57 
together with the three spinodals delineating the limits 
of metastability of the demixed fluid, mixed fluid and 
vapour (marked Si, S2 and 6*3). Also shown is the 'hid- 
den binodal' for coexistence between a vapour and a 
demixed liquid, calculated by neglecting the coupling of 
the concentration to the density. Clearly for this value of 
S, the hidden binodal and the liquid-vapour critical point 
(at which it terminates) both lie within the metastable 
region. For smaller S ^ 0.45, however, the metastable 
critical point moves outwith the limit of metastability 
and the hidden binodal is lost [ig. ||b]. Although the 
phases corresponding to the hidden binodal are not ob- 
servable at equilibrium, they are expected to influence 
the dynamical properties of the system. Wc return to 
this point in more detail in subsection [II B , 
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FIG. 3. (a) Liquid- vapour phase envelope in the p-T plane 
for 5 = 0.57. Also shown are the spinodals S\ , 5* 2 and 5*3 
and the 'hidden' binodal as described in subsections 
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(b) The phase envelope and spinodals for 5 = 0.45, by 



which point the metastable LV critical point has been lost. 



B. General Landau theory considerations 

We now turn to analyse the system from a more general 
point of view. Clearly we are dealing with a two order 
parameter problem, i.e. the density p = N/V and the 
number difference order parameter m — {Na ~ Nb)/N. 
In a symmetrical fluid, the Hamiltonian has to be invari- 
ant under sign reversal of m. The Landau expansion of 
the grand potential thus takes the general form 



F = a 



Kp- Po) 



(3.11) 



, m m B J , . 
+A— + — -jm^ip-po), 



where p is the chemical potential, and po is a reference 
density in the liquid-vapour coexistence region, chosen 
such that the cubic term oc {p — Po)^ vanishes. An expan- 
sion of this type has been discussed in a different context 
by Roux et aZ Pi; it applies generally to fluids with an 
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additional Ising like ordering tendency. In the case of 
Heisenberg type ordering, e.g., in a ferromagnetic fluid, 
the Landau expansion looks very similar with m simply 
replaced by the vector m, and many of the conclusions 
drawn below still hold. The liquid-vapour critical point 
of the mixed fluid is found at a = 1, and could the fluid 
be kept at fixed density po, it would demix or order at 
A = 0. 

Phase diagrams of such fluids are often discussed in 
terms of an 'interaction ratio' R, comparing the strength 
of the ordering or demixing tendency in the fluid with the 
overall attractive interactions between particles |0,^,|lO|. 
R corresponds to 1 — (5 in our model, to J/K in the 
Blume- Emery-Griffiths model ||3^, to the dipole moment 
m in reduced units in ||l^ and to 1/T° in Reference [0. 
It is generally found that increasing R drives the phase 
behaviour from the topology depicted in fig |l|a via 
towards i.e. a critical end point turns into a tri- 
critical point, which moves up in temperature, until the 
liquid-vapour critical point disappears in the region of 
coexistence between demixed and mixed liquid phases. 

Although the interaction ratio appears to be an influ- 
ential quantity, it is important to note that the basic fac- 
tor driving the phase behaviour is the c oupling between 
the two order parameters. In eq (3.11), it is described 
by the last term. At B~Q, no tricritical point can be 
expected, regardless of the interaction ratio. Therefore, 
in what follows we shall analyse the phase behaviour in 
terms of the coupling strength rather than the interac- 
tion ratio. The relation between the two will be discussed 
later. 

We will focus on the transition from a topology with 
a critical end point to one with a tricritical point, and 
assume that we are well below the liquid-vapour critical 
temperature, a < 0. It is the n convenient to rescale the 
Landau expansion eqn ( ^.llD such that 



Thus the A line shifts at the LV coexistence line. Techni- 



F = 



- 4 9 4 

— -\ + ^-fiQ + K{\-e)rn\ (3.12) 



where fi = /^/(a 



g = {p- po)/\/^, m 



and F is written in units of a . The parameter k = 
B/2^/—a then describes the effective coupling between 
the order parameter and the density. The parameter 
9 = A/ {—a) — 2k is temperature-like: 9 oc {T — Tcep), 
where Tcep is the temperature of the (stable or unsta- 
ble) critical demixing point at fixed density g — 1. The 
phase behaviour is found by minimising F with respect 
to m and g. 

At K = 0, the ordering behaviour and the liquid-vapour 
phase separation decouple and do not affect one another. 
One finds coexistence of a liquid {g — 1) and a vapour 
(^i = —1) at /i = 0., and these phase boundaries are 
crossed by the critical ordering (A) line at 6* = 0. If a 
small coupling k is turned on, the ordering temperature 
increases with the density: 



cally, the Landau expansion (3.12) predicts two critical 
end points, one on the liquid side at 6* = 0, and one on the 
vapour side a,t 9 = —4k. In all fluids studied so far, only 
the upper one has been seen. Situations with two CEPs 
are however encountered when a critical line intersects a 
liquid-solid coexistence region |0,|l^. 

Next we study the stability of the demixed liquid 
phase. The order parameter in the homogeneous demixed 
liquid takes the value = 9c{g) — 9. The determi- 
nant of the stability matrix there is given 

by 4m2((3g,2_ 

l)/2 — K^). Hence the demixed phase becomes unstable 
for g < gc, with the spinodal line 



gc = V(2k' + l)/3. 



(3.14) 



A tricritical point is found when the spinodal line inter- 
sects the critical line ( p.l3| ). This requires gc > 1, i.e. the 
coupling K has to be larger than a limiting value kq = 1. 
The tricritical point is then located at 



(gc) and gt = gc 



2 1 + 2k^ 



or M*=3 



(.2-1). 



(3.15) 
(3.16) 



We conclude that the coupling k between the order pa- 
rameter and the density determines the topology of the 
phase diagram. If k exceeds kq, the topology switches 
from one with a critical end point [fig. |l]a] to one with a 
tricritical point [fig. |l|b] . From a physical point of view, k 
basically refiects the correlations between order parame- 
ter and density. For example, the response of the average 
order parameter to a change of chemical potential in the 
demixed liquid phase is given by 



dm 
dp. 



(X (mg) — {m){g) cx k 



g~ gc 



(3.17) 



9c{g) = 2K{g-l). 



(3.13) 



where the exponent C is given by C = 1/2 in mean field 
theory. (Scaling arguments yield C — (7 + «)/2, where 7 
and a are the usual Ising critical exponents of the order 
parameter susceptibility and the specific heat). 

It is instructive to investigate the relationship between 
K and the interaction ratio R. When R increases, the 
critical end point temperature Tcep moves closer to 
the liquid- vapour critical temperature Tcq. Assuming 
a cx (Tco — T), this implies that the effective coupling 
K cx y/l/ (Tco — T) increases also, and diverges as Tcep 
approaches Tcq. As long as there is any coupling be- 
tween the density and the order parameter (i.e. B > 
in eq.( 3.1l| )), tuning R has the effect of tuning the cou- 
pling. Thus our arguments are supported by our own ex- 
plicit model calculations and by the results quoted earlier 

Next we discuss the implications for the LV phase 
boundary. The small coupling limit n < kq has been 
studied in detail in ref. [ psf . In chemical potential space, 
the critical end point induces a weak singularity in the 
first order liquid- vapour line |2q]. 
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fliT) = ^lreg{T) ~ U\t\' 



(3.18) 



where firegiT) is an analytical function of the tempera- 
ture T, t = {T — TcEp) /TcEP and U is & critical am- 
plitude. In the Landau theory framework, ^reg is simply 

given by ^Ireg = Areg = 0. 

The density of the vapour phase at coexistence behaves 
in a similar way 



Pg{T)^Pg^reg{T)~VgW 



(3.19) 



whereas the density of the liquid phase shows the marked 
hump indicated in fig. |^a, given by ]25t 



pi{T)^ pUreg{T)~Vl\t\'' 



(3.20) 



In the strong coupling limit, k > kqj one deals with the 
simple case in which two first order lines meet at a triple 
point. One thus expects a kink in p{T) and in Pg{T), 
and the density of pijT) jumps to the coexisting demixed 
liquid phase [cf. fig. ||c]. The two regimes meet at k = kq, 
where the critical end point turns into a tricritical end 
point (TEP) ^ (cf. fig|b). Within our Landau theory, 
we have calculated the phase boundaries in the vicinity 
of a tricritical end point. Since we expect that k varies 
with temperature, we consider a path k — kq + K9, where 
9 (X {T — Ttep) as defined above. Our results to leading 
order in 6 are summarised as follows: For the chemical 
potential of the liquid vapour line, we find 

/1 = (61 >0), (x^ ~{\9\/if/'^ [9<Q). (3.21) 

The rescaled densities of the vapour phase (spectator 
phase, Q-) and the liquid phase {q+) are given by 



Q- 
Q+ 



= 1 



= Tl (^ > 0) 

- 1/2(1^1/3)3/^ ( 
f (|0|/3)V2 ^ 



< 0). (3.22) 



Note that these results agree with the arguments pre- 
sented in 1^ if one inserts at = 1/2, the mean field 
value of the exponent a at the tricritical point. Mean 
field theory is expected to yield the correct tricritical be- 
haviour in our system, since the upper critical dimension 
of a tricritical point is d„ = 3. 

We close this subsection with a discussion of 'hidden' 
part s of the phase diagram. To this end, we return to 
eqn 3.11 and perform a general stability analysis. We 
consider the demixed phase at to = [A — B{p — po)]^^'^ 
and the mixed phase at to = 0. 

The conditions for stability of the demixed phase have 
already been discussed earlier. One finds a spinodal at 



p-po = ±Si with 5i = V(S2 - 2a)/6, (3.23) 



which is equivalent to eqn. 3.14 in non rescaled units. For 
\p — po\ < SI, the demixed phase becomes unstable with 



respect to phase separation into a demixed liquid and a 
vapour. 

In the present context, the stability of the mixed liquid 
phase is more interesting. The stability analysis yields 
two spinodals: At 



\p-Po\<±s; , 52*-y=«73, 



(3.24) 



the demixed liquid phase is unstable with respect to 
phase separation into a demixed liquid and its vapor, 
and at 



(3.25) 



it becomes unstable with respect to demixing. As long 
as ^ > and the coupling B is sufficiently small, there 
exists a region on the high density side of the fluid, p > 
POi where the mixed liquid phase can be metastable or 
even stable 

\/-a/3 < p- po < A/B 

In that binodal can be found at \p — po\ — y/~a, 

which may be stable or hidden in the metastable region. 

The spinodals and the hidden binodal are indicated in 
fig. 1^. Note that both the coefficients a and A depend 
roughly linearly on the temperature T. This explains 
the linear form of S2 in the density-temperature plane as 
opposed to the parabolic form of S3 . 

The hidden binodal disappears completely as soon as 
A < at a = 0, i.e., as soon as the spinodal 5*3 meets S2 
at the 'hidden critical point' or beyond (at p < po). This 
criterion is independent of the coupling B. The strength 
of coupling thus has no infiuence on the appearance of 
a hidden binodal. It does however affect the range of 
metastability of the hidden binodal, i.e., the temperature 
interval before it is lost by intersecting the spinodal 6*2. 

Of course, the concepts of spinodals and metastability 
only really make sense within a mean field treatment, and 
strictly speaking lose their physical meaning as soon as 
fiuctuations are taken into account. However, a discus- 
sion of the metastable and unstable regions is still useful 
in the context of the dynamical properties of the system. 
Consider, for instance, a binary fluid with interactions 
corresponding to the situation shown in fig. |^a, at density 
p = Pi), which is quenched from some high temperature 
into the coexistence region slightly below the 'metastable 
critical point'. One can then expect 'two-stage demixing'. 
In the first stage, the fluid will separate into domains of 
vapour and mixed liquid, and the separation of these do- 
mains will be accelerated by the driving force of gravi- 
tation. In the second stage, the liquid phase will slowly 
demix, and droplets of demixed liquid will additionally 
nucleate from the vapour phase. If the interactions of 
the fluid correspond to the situation depicted in flg. ||b, 
on the other hand, no intermediate stage will appear and 
the fluid will demix and phase separate simultaneously. 
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IV. MONTE-CARLO STUDIES 
A. Simulation details 

Many features of the simulation techniques employed 
in the present study have previously been detailed else- 
where [ 3^,^.0 . Accordingly, we confine the description 
of our methodology to its barest essentials, except where 
necessary to detail a new aspect. 

We assume our system to be contained in volume 
and to be thermodynamically open so that the total num- 
ber density and concentration can fluctuate. The associ- 
ated (grand canonical) partition function is: 



N 



EEn 

N=0{si} j=l 



dfi 



-0[^l>l{r,s}}+^iN] 



(4.1) 



Here, the species label Si — 1,-1 denotes respectively 
the two particle species A and B, N = Na + Nb is the 
total particle number, /3 = l/fc^T is the inverse tempera- 
ture and ^ is the chemical potential. The configurational 
energy density $ is given by 



(4.2) 



where the symmetri cal s qua re- w ell interparticle potential 
U is defined in eqs. and |2.2| . 

Grand canonical Monte Carlo (MC) simulations were 
performed using a standard Metropolis algorithm |^,|^ . 
The MC scheme comprises two types of operations: 



of simulations, we employed the histogram reweighting 
technique |39| . Use of this technique permits histograms 
obtained at one set of model parameters to be reweighted 
to yield estimates appropriate to another set of model pa- 
rameters. To enable simultaneous reweighting in all three 
fields fi, T, S, one must sample the conjugate observables 
{p,u,Ud), with p — N/V the number density, u = ^/V 
the configurational energy density, and Ud that part of u 
associated with interactions between dissimilar particle 
species. In addition to these variables, we have also ac- 
cumulated the quantity rh = {Na — Nb)/V = pm which 
gives a measure of the degree oi A — B ordering in the 
system. 

As mentioned in the introduction, standard GCE sim- 
ulations, deep within the LV coexistence region, are ham- 
pered by the large free energy barrier separating the two 
coexisting phases. This barrier leads to metastability ef- 
fects and prohibitively long correlation times. To circum- 
vent this difficulty, we have employed the multicanonical 
preweighting method | |2^ which encourages the simula- 
tion to sample the interfacial configurations of intrinsi- 
cally low probability. This is achieved by incorporating a 
suitably chosen weight function in the MC update prob- 
abilities. The weights are subsequently 'folded out' from 
the sampled distribution to yield the correct Boltzmann 
distributed quantities. Use of this method permits the 
direct measurement of the distribution of observables at 
first order phase transitions, even when these distribu- 
tions span many decades of probability. Details concern- 
ing the implementation of the techniques can be found 
in references [B2,n7|. 



1. Particle insertions and deletions. 

2. Particle identity transformations: A ^ B, B ^ A 

Since particle positions are sampled implicitly via the 
random particle transfer step, no additional translation 
algorithm is required. 

To simplify identification of particle interactions, we 
employed a linked list scheme |38| . This involves parti- 
tioning the periodic simulation space of volume into 

cubic cells, each of linear dimension the interaction 
range, i.e. L/l — 1.5. We chose to study two system 
sizes corresponding to I = 8 and / = 10, containing, at 
LV coexistence, approximate average particle numbers 
(N) — 550 and 1100 respectively. Equilibration periods 
of up to 2 X 10^ MCS were employed and sampling fre- 
quencies were 100 MCS for the / = 8 system to 150 MCS 
for the / = 10 system. Production runs amounted to 
2 X 10^ MCS for the Z = 8 and 5 x 10^ MCS for the I = 10 
system size. At coexistence, the average acceptance rate 
for particle transfers was approximately 10%, while for 
spin flip attempts the acceptance rate was approximately 
40%. 

In this work we wish to explore the parameter space 
spanned by the three variables (/x,T, (5). To accomplish 
this, without having to perform a very large number 



B. Method and results 

Using the multicanonical simulation scheme, we have 
obtained the density distribution p{p) for a number of 
states close to the LV coexistence curve, and for a number 
of choices of i5. We begin, by probing the regime of critical 
end point behaviour. 

On the basis of the mean field results of sec. Ill, CEP 
behaviour is expected to occur for large 5 < 1. For 6^1, 
the CEP will occur at very low temperatures relative to 
the LV critical point (i.e. T Tcq), but will move to 
higher temperatures as S is reduced. At some point as 6 
is reduced, the phase diagram is predicted to evolve into 
a triple point topology. Thus in seeking to observe CEP 
behaviour one should aim to set 6 large and to search at 
low temperatures. Unfortunately, very low temperatures 
are associated with high liquid densities at LV coexis- 
tence, and these inaccessible to our GCE scheme due to 
the prohibitively small particle transfer acceptance rate. 
In fact we find that the largest value of S for which the 
density fell within the accessible range {p ^ 0.7) was 
S = 0.72. Although this value of S is not as large as 
one might hope to attain, it nevertheless transpires that 
CEP behaviour occurs. This is demonstrated in fig. 
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which shows the hquid and vapour coexistence densities 
for this value of d obtained as the first moment of the re- 
spective peaks of p{p) for the I — 8 s yst em size The 
data were obtained by reweighting [p9| four histograms 
spanning temperatures in the range T — 0.99— 1.05, and 
coexistence was located using the equal peak-weight cri- 
terion for p{p) pO| , p4t . Also shown in the figure is the 
measured locus of the critical line. 



1.09 
1.07 
1.05 
1.03 
1.01 
0.99 
0.97 
0.95 

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 
P 

FIG. 4. The liquid-vapour coexistence curve in the p-T 
plane for 5 = 0.72, showing the vapour (V), mixed fluid (MF) 
and demixed fluid (DF) phases. The results were obtained 
from the measured peak positions of the coexistence density 
distributions for the / = 8 system size. Statistical errors do 
not exceed the symbol sizes. 

Clearly the data of fig. ^ display an anomaly in the liq- 
uid branch density close to the intersection point of the A 
line and the liquid branch {i.e. at the CEP). In the ther- 
modynamic limit, the liquid branch density is expected 
to exhibi t a cu sp-like singularity at the CEP, as given in 
equation 3.2C. In our finite-sized system, however, this 
critical singularity is smeared out and shifted, so that 









\ MF 










I? 









only a rounded depression in the coexistence envelope is 
visible Since these aspects of the CEP singulari- 

ties have recently been discussed in detail elsewhere , 
we shall not pursue them further here. Instead we shall 
proceed to consider what happens as d is made smaller 
still. 

Further reducing 6 continues to shift the CEP closer to 
the LV critical point. As one reaches 5 — 0.675, however, 
the phase diagram changes topology. We find that above 
a certain temperature the liquid peak inp^p) decomposes 
into two peaks. This is shown in fig. g for the / — 10 
system size at a temperature T = 1.044. Evident from 
this figure are two closely separated overlapping peaks, 
the presence of which signifies incipient triple point be- 
haviour. It follows that for this 6 and T, the system lies 
close to the tricritical end point which heralds entry into 
the triple point phase diagram topology [cf. fig. |^] . Ac- 
tually we believe the TEP lies close to 6 = 0.68 since this 
is the value at which we first observe the appearance of a 
shoulder in the liquid peak. In a sufficiently large system. 



this shoulder would presumably resolve itself into a dis- 
tinct peak. We have not, however, attempted to pinpoint 
the location of the TEP more precisely, as this would re- 
quire a full finite-size scaling analysis~a task beyond the 
scope of the present study. 




FIG. 5. The measured near-coexistence density distribu- 
tion p{p) for T = 1.044, S = 0.675 showing the three-peak 
structure discussed in the text. The distribution is normalised 
to unit integrated weight and statistical errors are comparable 
with the symbol sizes. 

Using histogram reweighting, we have monitored the 
temperature dependence of p{p) as 6 is reduced below 
the value at which the TEP occurs, fig. ||a shows a se- 
lection of density distributions for S = 0.665, for which a 
triple point (TP) occurs at some temperature Ttp < Tco, 
i.e. below the liquid- vapour critical temperature. The 
corresponding forms of p{rh) are shown in fig. At 
the triple point, a demixed liquid coexists with a mixed 
liquid and its vapour. For T > Ttp, there is phase coex- 
istence either between the mixed liquid and its vapour, or 
between the mixed and demixed liquid [cf. fig. The 
liquid-vapour coexistence terminates at the LV critical 
point, while the mixed-demixed liquid coexistence curve 
terminates at a tricritical point. From fig. ^ one sees 
that for 5 — 0.665, the tricritical point temperatures lies 
slightly below the LV critical point temperature, as evi- 
denced by the fact that on increasing T the liquid peaks 
merge before the liquid and vapour peaks do so. 
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FIG. 6. (a) Coexistence density distributions p{p) for 
5 — 0.665 at a selection of temperatures spanning the triple 
point temperature, (b) The corresponding form of p{rh) 
where ifi = mp = {Na — Nb)/V. Lines are guides to the 
eye. Statistical errors are comparable with the symbol sizes. 

The coexistence density distributions for S — 0.66 are 
shown in fig ^ for temperatures spanning the triple point 
temperature. For this S, the tricritical point is sufficiently 
well separated from the LV line that it is possible to 
distinguish the liquid-vapour and liquid-liquid branches 
by appropriately tuning the chemical potential. This is 
demonstrated in fig. ^ which shows p{p) for the two 
coexistence curves at T = 1.058. The different degree 
of order in the two liquid phases is clearly seen in the 
distribution p{m), shown in fig. One notices, how- 
ever, that both the density distributions show signs of 
the third phase. This reflects the closeness of the two 
coexistence curves at this S and T, as evidenced by the 
very small chemical potential difference. Under such con- 
ditions, finite-size smearing effects render it difficult to 
completely isolate two of the three phases. 
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FIG. 7. (a) Selected coexistence density distributions for 
S = 0.66, spanning the triple point, (b) The density distri- 
bution for T = 1.058 for two different values of the chemi- 
cal potential, (c) The corresponding forms of p(m), where 
rh — mp = {Na — Nb)/V. Lines are merely guides to the eye 
and statistical errors are comparable with the symbol sizes. 

Finally in this section, we consider the phase behaviour 
for S — 0.65. Coexistence forms of p{p) at selection of 
temperatures are shown fig ^. One observes that on in- 
creasing temperature, the low density vapour peak moves 
smoothly over to merge with the high density peak of 
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the ordered liquid. At no point is a three-peaked struc- 
ture visible. This scenario is consistent with the phase 
behaviour shown schematically in fig. |^d, in which the 
vapour and the demixcd liquid phases merge at a tricrit- 
ical point. 
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FIG. 8. Selected near-coexistence density distributions for 
S = 0.65 at a number of sub-tricritical temperatures. Statis- 
tical errors are comparable with the symbol sizes. 



V. DISCUSSION AND CONCLUSIONS 

In summary, we have used multicanonical Monte Carlo 
simulations and histogram reweighting techniques to 
study how the liquid-vapour phase behaviour of a sym- 
metrical binary mixture depends on S, the ratio of in- 
teraction strengths for dissimilar and similar particle 
species. For ^ ^ 1, the phase diagram exhibits a critical 
end point at temperatures well below the liquid-vapour 
critical point. Decreasing 6 shifts the critical end point 
closer to the liquid-vapour critical point. For S « 0.68, 
however, the critical end point becomes locally unsta- 
ble, and a triple point occurs in which vapour, a mixed 
liquid, and a demixed liquid all coexist. For tempera- 
tures above the triple point there is coexistence either 
between a high density demixed fluid and a moderate 
density mixed fluid, or between a mixed fluid and its 
vapour. Decreasing S still further pushes the triple point 
to higher temperature, until for 6 < 0.65 it eventually 
equals that of the isotropic liquid- vapour critical point. 
Thereafter, the mixed liquid phase is preempted by the 
demixed liquid phase and the liquid-vapour coexistence 
curve terminates in a tricritical point. 

Thus our simulation results confirm the qualitative pic- 
ture of phase diagram topology emer ging from mean field 
theory as set out in sections || and VU. These theories 
seem quite successful in capturing key features of the 
behaviour such as the existence of a coexistence curve 
anomaly at the CEP, the existence of the triple point 
regime, and the crossover to a purely tricritical regime 
Additionally, our Landau theory study of subsection 



the coupling of density and concentrations leads to the 
observed phase behaviour. 

In quantitative terms, however, the mean field theories 
are less reliable. Owing to the neglect of correlations, 
they predict neither the correct exponents for the coex- 
istence curve singularities at the CEP, nor the shape of 
the near critical LV coexistence curve. The values they 
yield for quantities such as the LV critical temperature 
are also at variance with simulation estimates: e.g. for 
S = 0.72 mean field calculations predict that the LV crit- 
ical temperature is T^q'^ = 1.172, while simulation gives 
^co™ = 1.06(1). In view of this, the apparently better 
agreement between the mean field and simulation esti- 



mates of the CEP for 5 = 0.72, 

rpsim 
^CEP 



T, 



mf 



CEP 



[II B 



provides useful physical insight into the manner in which 



1.002 and 

1.02(1), are presumable fortuitous. For 3d tri- 
critical behaviour, mean field theory is at least expected 
to yield the correct tricritical exponents, since the up- 
per critical dimension for such behaviour is d = 3 p5| . 
Although we have not attempted to probe the univer- 
sal aspects of the tricritical behaviour, our results show 
that estimates for the tricritical temperature are not re- 
produced by the simulations, at least close to the triple 
point regime. However, this may partly be the result of 
crossover effects associated with the relative proximity of 
the LV critical and tricritical points. 

The mean field estimates are also inaccurate regard- 
ing the sensitivity of the phase diagram topology to 
changes in 5. The calculations of section [II A| predict 
that the regime of triple point topology lies in the range 
0.605 < 5 < 0.708. In contrast, the simulation re- 
sults show this range to be considerably smaller, namely 
0.65 ^ (5 ^ 0.68. Indeed were it not for our use of his- 
togram reweighting to scan the phase behaviour as a 
function of (5, we might easily have missed this regime 
altogether. Thus it seems that more sophisticated liq- 
uid state theories are called for before the goal of accu- 
rately predicting the phase behaviour of simple binary 
fluid models is attained. Presumably any successful the- 
ory must be capable of dealing both with the critical and 
non-critical regimes of the phase diagram. Infact one 
such theory, the Heirarchical Reference Theory has 
recently been proposed. It would be interesting to see 
whether or not it accurately reproduces the phase be- 
haviour of the present model. 

With regard to further work on this model, one par- 
ticularly interesting project would be to investigate the 
predicted existence of the hidden binodal and associated 
metastable critical point. The occurence of metastable 
critical points was first discussed by Cahn |33|, and later 
found in lattice gas models by Hall and Stellpl . More re- 
cently it has suggested that they occur in colloidal fluids 
close to the freezing line [Q, in dipolar fluids and 
in models for water Since the present model offers 
a computationally tractable system, it might usefully be 
employed as a test-bed for studying the generic features 
of the metastable critical point. This could feasibly be 
achieved by quenching the system from high temperature 
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into the unstable regime just bel ow th e metastable criti- 
cal point. As described in section III E| , this should result 
in a two-stage demixing process in which the metastable 
mixed liquid phase appears for a transitory period before 
eventually demixing at later times. 

Additional interesting work would be to look at the 
equilibrium phase behaviour of the symmetrical mixture 
as a function of S < 0. Landau theory Q predicts that 
as S is made increasingly negative, the tricritical point 
transforms first into a double critical end point, before a 
critical end point emerges on the vapour side of the LV 
coexistence envelope. It would certainly be worthwhile 
to assess whether or not this scenario is correct. 
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